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Abstract 


The cumulative degree distributions of transport networks, such as air transportation net- 
works and respiratory neuronal networks, follow power laws. The significance of power 
laws with respect to other network performance measures, such as throughput and syn- 
chronization, remains an open question. Evolving methods for the analysis and design of 
air transportation networks must be able to address network performance in the face of 
increasing demands and the need to contain and control local network disturbances, such as 
congestion. Toward this end, we investigate functional relationships that govern the perfor- 
mance of transport networks; for example, the links between the first nontrivial eigenvalue, 
A 2 , of a network’s Laplacian matrix — a quantitative measure of network synchronizability — 
and other global network parameters. In particular, among networks with a fixed degree 
distribution and fixed network assortativity (a measure of a network’s preference to attach 
nodes based on a similarity or difference), those with small A 2 are shown to be poor syn- 
chronizers, to have much longer shortest paths and to have greater clustering in comparison 
to those with large A 2 . A simulation of a respiratory network adds data to our investigation. 
This study is a beginning step in developing metrics and design variables for the analysis 
and active design of air transport networks. 


1 Introduction 

The U.S. and world-wide air transport networks are scale-free 1 , i.e., their degree distributions 
follow power laws [9, 10]. A large body of recent work has been published on scale-free 
networks including the popular books Linked by Barabasi [4] and Six Degrees by Watts [21]. 
A number of excellent review articles, including Newman [17], Strogatz [19], Albert and 
Barabasi [1] and Dorogovtsev and Mendes [6], contain hundreds of references. 

Arguably, many scale-free networks occurring in natural and technological realms have 
never been actively designed in the traditional sense: identify design variables, objectives, 
and constraints and follow a prescriptive algorithm to obtain a design that satisfies con- 
straints and is “optimal” with respect to the given objectives. Rather, these networks have 
evolved in response to demands, in accordance with some natural or technological rules. 
Thus many of the investigations have been of an analytical nature, i.e., given a particular 
natural or technological network, its characteristics are studied. Our ultimate interest is 
to take a step from analysis to active design and our motivation comes from air transport 
systems. We emphasize that complex networks will likely never be completely amenable to 
traditional design methods, given the intrinsic lack of predictive modeling akin to that of 
physical artifacts (e.g., airplanes, automobiles). However, we conjecture that some measure 
of active design is still possible with the identification of appropriate design variables and 
metrics. Specifically, we are looking for an appropriate functional relation between global 
metrics (e.g., throughput, delays, capacity, synchronizability) and locally controllable struc- 
ture (e.g., connectivity, degree). 

Air transport networks are our ultimate domain of interest. Are scale-free networks 
desirable for air transport? Given a fixed degree distribution, how should the network links 
be configured to achieve optimal performance for relevant metrics? Can salient features of 
scale-free air transport networks critical to network performance be identified? 

With the advent of deregulation of the U.S. air transport system in 1978, airlines began 
to organize their operations with a hub-and-spoke approach. Two natural outcomes have 
been an increase in flight frequency and an increase in the variance of flight times. Recently, 

1 We use the term “scale-free” here to denote networks whose degree distributions follow a power law. 
Due to a wide range of properties possessed by networks of similar degree distributions, there is an ongoing 
discussion of the meaning of “scale- free” [13]. 
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competitors to the hub-and-spoke model have garnered attention with the use of point- 
to-point flight schedules. With the delays experienced by travelers at hub airports these 
direct flights are an attractive alternative. How will these point-to-point airlines alter the 
air transport network structure? Can we provide any guidance to local air routing decisions 
with global air transport network performance measures in mind? To begin to address these 
questions we examine the effects of a network metric for synchronization on transport route 
structure by holding the degree distribution and a scale-free/scale-rich metric constant. 
We also report on a simulation of a respiratory neuronal network used as an additional 
testbed for investigating the synchrony metrics. It is a first step in the investigation aimed 
at deriving the functional dependences among various local network properties and the 
aggregated metrics of interest to participants in the transportation system in an effort to 
eventually arrive at active transport network design and optimization algorithms. 


2 Background 

In this section we briefly review some of the network attributes salient to our investigation 
into the network functional relationships. Early network growth models were based on 
preferential attachment. A variety of authors have developed extensions and improvements 
to the early models. All of these mechanisms build networks sequentially, one node at a 
time. Barabasi’s [4] original approach selected the end nodes for the edges associated with 
the new node based solely on the current network’s degree distribution. Subsequent efforts 
have altered the end node selection method to control other network features. For example, 
Wang et al. [20] develop a growth model in which assortativity is tunable while Schank and 
Wagner [18] and Holme and Kim [11] grow networks with tunable clustering coefficients. 

The network attribute under study here is synchronicity. One notion of synchronicity has 
to do with the network’s tendency to synchronize over time, given a specific static structure 
of the network. In particular, we are not yet considering explicit traffic flows through the air 
transport network. Instead, we are investigating how the node (e.g., airport) connectivity 
may influence the traffic flow. 

Of importance here is the tunability of a given network with respect to synchronicity. 
We realize tunability as network optimization. Before proceeding further we provide a 
definition of network synchronization for a discrete complex system. Given a connected 
network, denote the state of a node i at time t by Xi(t). How do the states of the nodes 
change over time? Clearly if nodes do not rely on any information generated by adjacent 
nodes then there is no opportunity of synchronization. Atay et al. [2, 3] assume that all 
nodes are identical and conform to the following generic discrete time equation to determine 
their next state: 


Xi(t + 1 ) = f(Xi(t)) + K 


X /(**(*)) - 

’ j\(hj)£Edges 


(1) 


where k. known as the coupling constant , is a scalar describing the extent to which neighbors 
effect the state of a node; /(•) is any differentiable function mapping some finite interval 
to itself. The function /(•) describes the behavior of a node in the absence of any outside 
influence. We say that a network synchronizes for a given initial condition if for all i,j 


lim | Xi(t) — Xj(t) | = 0. (2) 

t— >oo 

Note that if k = 0, then the equation becomes Xi(t + 1) = f{xi(t)). 

One attribute that correlates with the network’s capacity to synchronize is the first 
nontrivial eigenvalue, A 2 , of the Laplacian matrix associated with the network structure 
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(more about A 2 in the next section). Here A 2 acts as a measure of the range of k over which 
the network will synchronize. 

Clearly, air transport networks will not meet some of the assumptions. The nodes 
(airports) are not all identical and we are more interested in the transient (say, 24 hour 
period) behavior than in what happens as t -> 00 . Nonetheless, in the next section we will 
see that A 2 will provide useful information with regard to network structure and synchrony. 

For some complex systems synchronization is an essential feature. For example, in 
mammals a small group of neurons (roughly 200) is responsible for generating a regular 
rhythmic output to motor cells that initiate a breath. (We explore this example further in 
Section IV.) Without synchronization of the neuronal output, breathing would be ragged or 
not occur at all. However, in the defined sense, synchronization is an undesirable attribute 
for air transport networks. Think of the airports as neurons in our mammalian respiratory 
example. Inhaling means all planes land at all airports simultaneously. Exhaling means 
they depart together. The result is congestion. Thus, for the given definition, one would 
like an air transport network design to minimize synchronization. 

Another network metric, developed in [13], serves to classify the ways in which networks 
with a given degree distribution may be constructed. The metric is s(G), where G stands 
for “graph” (networks and graphs are interchangeable). To determine s(G), compute the 
product of the degrees of the end nodes for each edge or link, sum them up for all edges, 
and divide by s max , where s max is the maximum value of the sum taken over all possible 
connected graphs for a fixed degree distribution. That is, 

^(CO — ( d^ x dj)/s max , (3) 

( i,j)EEdges ) 

where di is the degree of node i. The value of s max provides a way to scale the sum of 
the product of the degrees for each edge. In [13], for a fixed degree distribution, graph 
realizations with large values of s(G) are termed scale-free and graph realizations with small 
values of s(G) are termed scale-rich. Consequently, scale- free graphs or networks are those 
in which high degree nodes are more likely to be adjacent to other high degree nodes while 
scale-rich graphs are those for which high degree nodes are more likely to be adjacent to 
low degree nodes. As shown in [13], s(G) and assortativity are equivalent but are scaled in 
different ways. Typically the scaling for the assortativity measure leads to a tighter range 
of values for a given degree distribution. For further information about s(G) and how it 
is used to distinguish between networks for the Internet at the router level the interested 
reader is referred to [13] and [14]. 



Internet 

Air Transport 

product 

constraint 

competitors 

links 

distributors 

packets 

bandwidth 

ISPs 

hardwired 

routers 

planes (loaded) 
airport capacity 
airlines 
FA A/ Air lines 
airports 


Table 1. Analogy between Internet router and air transport networks. 


The problems faced by designers of air transport networks share some aspects with the 
design of an Internet router network. Many authors have contributed to investigations of 
how a router network is constructed. Two references in this field, [13] and [14], contain 
ideas central to our consideration of the design of air transport networks. At one level of 
resolution, Table 1 points out the analogies between these two network design problems. 
With regard to bandwidth, the Internet router designer must weigh the trade-offs between 
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Figure 1. Geometric graph with A 2 = 0.238 and A 2 = 0.925 


many low bandwidth connections and fewer high bandwidth connections. These trade-offs 
are akin to choosing between a few hub airports in a hub-and-spoke system and choosing 
lower frequency airports that might arise in a direct route system. Of course, there are many 
differences as well. The variation in the size of the packets for the Internet is not nearly 
as great as the number of passengers on planes of different sizes. In addition, although the 
FAA 2 clearly defines the routes allowed between airports, the links are as not hard-wired as 
they are in the Internet model. Still there is much to be learned from the research efforts 
on the design of effective Internet router networks. 


3 Network Experiments with A 2 

In this section the interplay between A 2 , a gross measure of network synchrony, and s(G), a 
network measure similar to assortativity is examined. Networks of two types — preferential 
attachment and geometric — are the testbeds. Each network has 100 nodes and is a simple 
undirected network (no self-loops, no multiple edges). An adjacency matrix, A, is con- 
structed. G(V. E) denotes a graph (network) with vertex (node) set V and edge (link) set 
E. For each network type the degree distribution and s(G) values are fixed. Finally, a tabu 
search heuristic rewires a given network so as to minimize (or maximize) the value of A 2 . 

Preferential attachment graphs are generated following the approach given in [4]. The 
network is grown by adding nodes and edges. For each node added, m edges are added 
preferentially, based on the current degree distribution. Geometric graphs are generated by 
randomly selecting 100 points (r, 6) with values of r £ [0, 1] and values of 9 € (0, 360]. Edges 
exist between pairs of points if the Euclidean distance is less than a specified threshold (in 
our experiments thresholds between 0.17 and 0.25 were used). If the resulting graph is 
connected, it is kept; otherwise it is rejected and the process begins again. A variety 
of network performance measures are available. These include network diameter, average 
degree, assortativity, clustering coefficient, synchrony and s(G). 

Our measure of synchrony relies on computing an eigenvalue of the Laplacian matrix 
associated with a given network structure. For our networks or graphs the Laplacian matrix 
is a symmetric matrix L = D — A, where D is a diagonal matrix with the degree of each 
node located along the main diagonal. A is the adjacency matrix for the graph. The 
second eigenvalue of L measures algebraic connectivity [7]. Intuitively, graphs with small 
A 2 are easier to “pull apart”. In particular, if A 2 = 0, then the graph is disconnected. 
Many authors [2,5,22] have convincingly used A 2 as a global measure of how likely a graph 
is to synchronize. That is, given an arbitrary flow of entities, the graph is less likely to 
synchronize if A 2 is small. The two graphs in Figure 1 have identical degree distributions, 
but the graph on the left is more weakly connected (e.g., the removal a single edge can 
disconnect the graph). The identification of structural differences between large graphs 

2 Federal Aviation Administration 
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with varying values of A 2 is studied in the following set of numerical experiments. 

In the remainder of the section, we describe numerical experiments in optimizing two 
types of networks for A 2 : geometric graphs with 100 nodes and preferential attachment 
graphs with m = 2. 

The graphs plotted in Figures 2 and 3 were constructed by first generating a random 
instance of the particular graph class — geometric in Figure 2 or preferential attachment 
in Figure 3. Next, a simple tabu search [8] heuristic was called to minimize or maximize 
A 2 while keeping the degree distribution and s(G) fixed. Allowable moves (re-wirings) are 
pair-wise edge interchanges that preserve the degree distribution and s(G). Briefly, the 
tabu search checks to see if the move is acceptable, that is, if the move is improving and 
not tabu, or improving and tabu but leads to the best observed value of A 2 (aspiration 
criterion). Note that these moves are precisely the moves allowed in a random rewiring 
scheme without checking for the preservation of s(G). The interested reader is referred to 
Glover and Laguna [8] for further information on tabu search. 

Figures 2-4 display networks with respect to the reciprocal of the eccentricity of each 
node u. The eccentricity of u is its maximum (shortest path) distance. The graphs are 
generated by socnetv 3 . The goal of the plots is to uncover any qualitative differences 
between the graphs with small and large values of the second eigenvalue of the Laplacian. 4 
For each pair of plots in a figure, the number of nodes (100), the degree distribution, and 
the value of s(G) is fixed. 

Nodes with equal eccentricity values are plotted on the same (dashed line) circles. The 
circles with larger radii have larger eccentricity. Consequently, nodes near the center have 
shorter longest paths. The paired plots exhibit large qualitative differences in the eccen- 
tricity patterns. The same pairs of plots for other available measures in socnetv were also 
constructed. Although small variations were seen, none of the other paired plots exhibited 
significant differences. 

Qualitatively, when A 2 is small, the patterns are less organized, the eccentricity plots in 
Figures 2a and 3a are more dispersed and consist of many rings of constant eccentricity. The 
eccentricity plots with larger A 2 are more organized, with few rings of constant eccentricity. 
Specifically, the plots with small A 2 have 11 (Figure 2a) and 10 (Figure 3a) rings. For 
larger A 2 there are five (Figure 2b) and four (Figure 3b) rings, respectively. The ranges of 
eccentricity values for the small A 2 plots are dominated by the ranges for the large A 2 plots. 
For example, the range of eccentricity values for the geometric graph with small A 2 (Figure 
2a) is [26,42] and for the geometric graph with large A 2 (Figure 2b) it is [4,8]. Thus, the 
patterns of eccentricities in Figures 2b-3b are non-overlapping and dominate those given in 
Figures 2a-3a. 

The diameter of the graph in Figure 2a is 42 while the graph diameter in Figure 2b is 
8. The diameters are 19 and 6 for graphs in Figure 3. Notice that this is also true for the 
simple graphs in Figure 1. The graph on the left has a larger diameter than the one on 
the right. For graphs with a fixed degree distribution and a fixed value of s(G), this result 
appears to hold in general. We know of no theorem that proves this result but numerous 
computational tests support this claim so far. Moreover, the inverse relationship between 
A 2 and the eccentricity does not hold if s(G) is allowed to vary. Figures 4a and 4b provide 
an example. The value of A 2 = 0.935 in Figure 4a is larger than A 2 = 0.440 in Figure 4b. 
Yet the range of eccentricity values for Figure 4a, [8, 16], is larger and does not overlap with 
the range for Figure 4b of [3,6]. 

In Figures 5 and 6, we investigate the relation between s(G) and the clustering coefficient, 
c(G) [17], of a network or graph. These figures display results for 5000 geometric graphs. 
In Figure 5, the 5000 geometric graphs are generated randomly. Edges connect nodes if the 

3 The source code and documentation can be found at http://socnetv.sourceforge.net/. 

4 We leave it to the reader to become acquainted the variety of measures and display features in socnetv. 
For the purposes of this exposition, we are interested only in the qualitative differences between the plots. 
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(a) A 2 = 0.009, e(V) = (26,42), c(G) = 0.426 (b) A 2 = 0.314, e(V) = (4,8), c(G) = 0.297 

Figure 2. Geometric graphs: 100 nodes, s(G) = 0.971, fixed degree distribution 



(a) A 2 = 0.006, e(V) = (10, 19), c(G) = 0.210 (b) A 2 = 0.365, e{V) = (3, 6), c(G) = 0.101 

Figure 3. Preferential attachment: 100 nodes, s(G) = 0.716, fixed degree distribution 


Euclidean distance between these nodes is less than 0.23. If the resulting graph is connected, 
it is kept; otherwise it is rejected and the process begins again. Figure 5b displays the values 
of s(G) versus c(G). There is no apparent correlation. Figure 6a records the values of s(G) 
when a given 100-node geometric graph is randomly rewired 5000 times. Figure 6b illustrates 
the inverse relationship between s(G) and c(G) when the first 364 re-wirings are excluded. 

In addition to the inverse relationship between \ 2 and the eccentricity, the clustering 
coefficient varies inversely with X 2 in Figures 2 and 3. Notice that when s(G) is not held 
constant, as in Figure 4, this relationship does not hold. A similar trend between s(G) 
and c(G) is observed in Figure 6. Figure 6a displays the value of s(G) versus 5000 random 
re-wirings of a given geometric graph. We note that the same moves as those used in our 
tabu search to optimize A 2 are used for the random re-wirings. Figure 6b plots s(G) versus 
c(G) for the last 4634 random re- wirings and s(G) varies inversely with c(G). Notice that in 
Figure 6a, the first 600 or so random re- wirings decrease s(G) almost monotonically before 
settling into an oscillating pattern of increases and decreases in the range of (0.79,0.83). 
This is not the case, however, when geometric graphs are generated at random (no rewiring), 
as in Figure 5. Here no correlation is exhibited between c(G) and s(G). 

As we have noted earlier, Figures 4a-b illustrate that s(G), or some yet unknown network 
metric, appears to exert significant influence on the eccentricity pattern. In Figure 2b, 
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(a) A 2 = 0.935, s(G) = 0.797, e{V) = (8,16), 
c(G) = 0.306 



(b) A 2 = 0.440, a(G) = 0.677, e(V) = (3,6), 
c(G) = 0.126 


Figure 4. Geometric graphs: 100 nodes, fixed degree distribution, varying s(G) 




(a) s(G) 


(b) s(G) vs. c(G) 


Figure 5. 5000 Random geometric graphs 


A 2 = 0.314, s(G') = 0.971, and the eccentricity range is [4, 8]. This compares favorably with 
the results in Figure 4b, where A 2 = 0.440 is larger and the eccentricity range of [3,6] is 
smaller with a smaller minimum node eccentricity. However, the comparisons with Figure 
4a are not consistent. The value for A 2 in Figure 4a is larger but, unexpectedly, the max 
and min values for the eccentricity range are much larger than those in Figure 2b. In 
addition, the previously observed pattern of c(G ) varying inversely with A 2 no longer holds. 
For example, A 2 decreases from 0.935 in Figure 4a to 0.440 in Figure 4b. Similar decreases 
observed in Figures 2 and 3 led to a doubling of c(G). But here c(G) decreases by more than 
a half. It is unclear if the role of s(G) explains the lack of consistency. Figure 5 provides 
one possible explanation. Here geometric graphs are generated at random (no rewiring) and 
there is no correlation between s(G) and c(G). In Figure 6, an inverse correlation exists 
but here the graphs are constructed by successively rewiring a single geometric graph at 
random. 
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Geometric Graph: 100 nodes, 5000 random rewirings 



Geometric Graph: 100 nodes, 4634 rewirings 
0.835 1 1 1 1 



0.795 

0 24 0 26 0 28 0 3 0.32 0 34 

Clustering Coefficient 


(a) s{G) 


(b) s(G) vs. c(G) 


Figure 6. 5000 Re- wirings of a Geometric Graphs 


4 Respiratory Network Simulation 

In the previous section we have seen that when the degree distribution and s(G) are fixed, 
there is a predictable difference in the shortest path distribution (eccentricity measures) and 
the clustering coefficients. Ideally the next step would be to simulate the performance of 
networks presented in the previous section and investigate their performance as air transport 
network systems. We are currently proceeding in this direction using a model previously de- 
veloped to simulate the complete U.S. air transport system. For now, we provide simulation 
results for a different system — a respiratory neural network — for which a simpler simulation 
model was readily available. In this model, synchronization (rhythmic breathing) is desired. 

Although synchronization is undesirable for air transport networks, there are systems 
for which it is an essential feature. In mammals, a small group of neurons is responsible 
for generating a regular rhythmic output to motor cells that initiate a breath. The network 
structure of these neurons allows them to synchronize without any external influence and 
produce regular bursts that lead to breaths. In [12], two geometric networks, one with a 
value of A 2 = 0.025 and a second with a value of A 2 = 0.974 were tested in a simulation model 
[16] of this neuronal network. The rhythmic output from the the network with A 2 = 0.025 
was ragged with fuzzy bursts, while outputs from the network with A 2 = 0.974 was sharp 
with clear, regular bursts (Figure 7). 

In mammals, a small group of neurons in the brain stem, called the pre-Botzinger com- 
plex, is responsible for generating a regular rhythmic output to motor cells that initiate 
a breath. Disconnected, these neurons are unable to provide sufficient output to activate 
the motor neurons, but their interconnected network structure allows them to synchronize 
without any external influence and produce regular bursts. Using a detailed simulation due 
by Hayes [16], we were able to experiment with how different network topologies control 
the effectiveness of the pre-Botzinger complex. We began by testing two geometric graphs 
with extreme values of A 2 . The results of the two simulations, depicted in Figure 7, provide 
compelling evidence for the utility of A 2 as a predictor of synchronization. It is easy to see 
that the network with higher A 2 synchronizes more strongly than the other network. The 
second set of simulations investigated two preferential attachment networks. The raster 
plots in Figure 8 are nearly indistinguishable. The results of the simulation are further 
analyzed via an autocorrelation analysis (Figure 8). Analysis (as in [16]) uncovers better 
synchronization in the network with the higher value of A 2 . The results in Figure 9 confirm 
that, although the difference is undetectable at a first glance (Figure 8), higher A 2 results in 


8 


Neuron spikes for X 2 = 0.97452 


Neurons spikes for X 2 = 0.02491 3 




Figure 7. Raster plots of neuron output for two networks with disparate A2 values. A 
point at (, x , y) indicates neuron x is spiking at time y. The higher A2 network displays much 
stronger synchronization among all nodes as predicted, as well as a quicker breath frequency. 


Neurons spikes for X 2 = 0.055883 


Neurons spikes for X 2 = 0.20587 




Figure 8. Indistinguishable raster plots of simulated neuron output for two sample networks 
with differing A2 values. 


a better synchronization. Autocorrelation indicates the largest difference during the refrac- 
tory (non-spiking) period: the two graphs exhibit similar behaviors during spikes, but not 
between spikes. These experiments provide further evidence that A2 can be used to identify 
graphs (networks) that are not likely to synchronize. 


5 Conclusions and Discussion 

Given the present state of air transportation networks, there is some urgency in developing 
active and rigorous design methodologies. Our goal is to develop a systematic way to 
design for some salient aspects of air transport networks. In particular, network structure, 
both static (node location) and dynamic (air route scheduling) has a direct effect on the 
functioning of the traffic in the network; we are now concerned with the effect of the static 
and dynamic network structure on the performance. 

How should the design process proceed? Design involves being able to manipulate vari- 
ables so as to optimize objectives subject to constraints. As a step in this direction we have 
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Figure 9. An autocorrelation plot of pre-Botzinger complex synchronization on two networks 
with the same degree distribution, but with differing A 2 values. The autocorrelation analysis 
shows that the higher A 2 network displays better synchronization. 


demonstrated that for a fixed degree distribution and fixed s(G) value (and, consequently, a 
fixed assortativity) , optimizing for A 2 yields networks with distinct eccentricity patterns. We 
have demonstrated the ability to construct networks with locally optimal A 2 and observed a 
correlation with global network attributes, such as clustering, eccentricity and synchroniz- 
ability. These results are further supported by a simulation analysis of another transport 
system - a respiratory neuronal network. This simulation supports our conjecture that large 
differences in A 2 result in observable differences in the burst activity: good synchronization 
for large A 2 and poor synchronization for small A 2 . It remains to simulate air transport net- 
works with small and large values of A 2 (more than likely with identical degree distributions 
and s(G) values). In addition to validating the static results for networks (as in Figures 2 
and 3) the simulation will also measure quantities of interest to the FAA that are currently 
not amenable to rigorous optimization. For example, the simulation will measure congestion 
effects in air traffic sectors. Finally, one of our major tasks is to derive maps between the 
network metrics we can control and airspace simulations and, ultimately, FAA metrics. We 
conjecture that deriving the maps will enable active design for a number of objectives and 
constraints. 
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